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The hydrodynamics of granular gases of viscoelastic particles, whose collision is described by 
an impact-velocity dependent coefficient of restitution, is developed using a modified Chapman- 
Enskog approach. We derive the hydrodynamic equations and the according transport coefficients 
with the assumption that the shape of the velocity distribution function follows adiabatically the 
decaying temperature. We show numerically that this approximation is justified up to intermediate 
dissipation. The transport coefficients and the coefficient of cooling are expressed in terms of 
the elastic and dissipative parameters of the particle material and by the gas parameters. The 
dependence of these coefficients on temperature differs qualitatively from that obtained with the 
simplifying assumption of a constant coefficient of restitution which was used in previous studies. 
The approach formulated for gases of viscoelastic particles may be applied also for other impact- 
velocity dependencies of the restitution coefficient. 
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I. INTRODUCTION 

Granular systems composed of a large number of dis- 
sipatively interacting particles behave in many respects 
as a continuous medium and may be, in principle, de- 
scribed by a set of hydrodynamic equations with appro- 
priate boundary conditions. Although this approach is 
successively used in various fields of engineering and soil 
mechanics (e.g. (3,0) a first-principle theory for dense 
granular media is still lacking. Hydrodynamics may be 
also applied to much simpler systems, such as rarefied 
granular gases. It has been used to describe many differ- 
ent processes, e^. rapid granular flows, structure forma- 
tion, etc. (see [Sj for an overview). For these systems the 
hydrodynamic equations are not postulated but derived 
from the Boltzmann equation. The corresponding trans- 
port coefficients are not phenomenological constants, in- 
stead they are obtained by regular methods, such as the 
Grad method 4] or the Chapman-Enskog method 
In most of the studies, which address the derivation of 
hydrodynamic equations and kinetic coefficients, it was 
assumed that the coefficient of restitution e is a material 
constant, e.g. |Ellllll|l3|lllIllI3. 

This assumption simplifies the analysis enormously, 
however, it is neither in agreement with experimental 
observations (e.g. 1141. llSl llq ) nor with basic mechanics 
of particle collisions The coefficient of restitution 

depends on the impact velocity g and tends to unity for 
very small g. As a consequence, particles behave more 
and more elastically as the average velocity of the grains 
decreases. The simplest collision model which accounts 
for dissipative material deformation, is the model of vis- 
coelastic particles. It is assumed in this model that the 
elastic stress in the bulk of the particle material depends 
linearly on the strain, whereas the dissipative stress de- 
pends linearly on the strain rate. 

Based on the fundamental work by Hertz the in- 
teraction force between colliding viscoelastic spheres has 



been derived [l3- Using this generalized Hertz law the 
coefficient of restitution of viscoelastic particles can be 
given as a function of the impact velocity and material 
parameters (20^: 



e = l~ CiAa'^^g'^' + \clA 



with 



g = |e • [vi - V2)\ = |e • vi2\ . 



(1) 



(2) 



The unit vector e — ri2 Iryi specifies the collision geome- 
try, i.e., the relative position r\2 = r\—r2 of the particles 
at the collision instant. Their pre-coUision velocities are 
given by v\ and V2- The elastic constant 



a = - 
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,cff 



depends on the effective mass and radius m" 
17111712/ {mi+m2), = i?ii?2/(i?i+-R2) of the colliding 
spheres, on the Young modulus Y, and on the Poisson 
ratio ly of the particle material. The dissipative coeffi- 
cient A is a function of dissipative and elastic constants 
(see ll9J for details). Finally, Ci is a numerical constant 



Ci = 



r(3/5) 



21/552/5 r (21/10) 



1.15344. 



(4) 



Equation ^ describes pure viscoelastic interaction. The 
assumption of viscoelastic deformation is justified if the 
impact velocity is not too large to avoid plastic defor- 
mation of the particles and not too small to neglect sur- 
face effects such as adhesion, van der Waals forces etc. 
We also assume that the rotational degrees of freedom of 
particles may be neglected and consider a granular gas of 
identical particles in the absence of external forces. Then 



2 



the coefficient of restitution gives the velocities of parti- 
cles after a collision w/, V2 in terms of their values before 
the collision: 



II. VELOCITY DISTRIBUTION AND 
TEMPERATURE IN THE HOMOGENEOUS 
COOLING STATE 



^1/2 = ^1/2 T 
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A. Evolution equations for temperature and for 
the second Sonine coefficient 



The impact velocity dependence of the coefficient of 
restitution implies serious consequences for the granu- 
lar gas dynamics: For the simplified case e = const, the 
form of the velocity distribution function f{v, t) is char- 
acterized by a time- independent scaled function /(c). It 
depends only on the scaled velocity c = v/vt, where 
vrii) = \/'^T{t)/m is the thermal velocity and T{t) is 
the granular temperature |22, HI] . The distribution 
function depends on time only via the time dependence 
of the granular temperature, its shape is time indepen- 
dent. The small deviations of the velocity distribution 
function from the Maxwell distribution are determined 
by the time-independent coefficient of restitution. For 
granular gases of viscoelastic particles, however, the ef- 
fective value of £ changes with time along with the ther- 
mal velocity vxit), which gives a typical velocity of the 
particles. Therefore, the shape of the velocity distribu- 
tion function evolves in a rather complicated way p3 | . 

The time dependence of the velocity distribution func- 
tion has to be taken into account when the hydrody- 
namic equations and the transport coefficients are de- 
rived. Since the simple scaling is violated, the stan- 
dard methods of kinetic theory of gases, developed for 
e = const, (e.g. 0,0,113) must be revised. 

As it follows from our analysis the transport coeffi- 
cients for gases of viscoelastic particles depend on tem- 
perature and time rather differently as compared with the 
case e = const. Correspondingly, the behavior of gases 
of viscoelastic particles differs qualitatively from that 
of gases of particles with the simplified collision model 
e = const. 

We wish to remark that the impact velocity depen- 
dence of £ has been already taken into account in Ref. 
|26l | for the hydrodynamic description of granular shear 
flow. In this study an empirical expression of the coef- 
ficient of restitution was applied and a Maxwellian ve- 
locity distribution was assumed. Moreover, the authors 
have used the standard Chapman-Enskog method with- 
out the modifications required for gases of dissipatively 
colliding particles. These modifications for the case of 
£ = const, have been extensively elaborated in [Tol |. 

The aim of the present study is to develop a continuum 
description of granular gases of viscoelastic particles. We 
derive the hydrodynamic equations along with the trans- 
port coefficients and the coefficient of cooling. In the rest 
of the paper, in Sec. II-V, we discuss in detail the most 
important case of the three-dimensional gases, while in 
Sec. VI we present the results for the two-dimensional 
systems, which are frequently addressed in molecular dy- 
namics studies. 



The Boltzmann equation for a granular gas of vis- 
coelastic particles in the homogeneous cooling state reads 
I2II23 



d_ 



f {vi,t) = cr'^92{a-) J dv2 J deQ {-V12 ■ e) \vi2 ■ e\ 

X [Xf{vl',t)f{4',t)-f{m,t)f{v2,t)] 
^ 92{<j)I{fJ) , (6) 



where a — 2R is the particle diameter. The velocities 
vl' and V2 denote the pre-coUision velocities of the in- 
verse collision, which leads to the after-collision veloc- 
ities vi and V2- The factor \vi2 ■ e \ characterizes the 
length of the collision cylinder of cross-section cr^ and 
the Heaviside step-function Q (— W12 • e) assures that only 
approaching particles collide. The contact value of pair 
correlation function 92 (cr) accounts for the increased col- 
lision frequency due to excluded volume effects. Finally, 
the factor x in the gain term accounts for the Jacobian of 
the transformation {vl',V2) {vi,V2) and for the ratio 
of the lengths of the collision cylinders jw/j ' 6*1 / l'S'12 • e | 
for the direct and the inverse collision. For the case of 
spheres colliding with a constant coefficient of restitution 
X — l/£^, while for viscoelastic spheres, with e = e{g) 
given by Eq. (|TJ, it reads p3.l27j 



X = 1 + ^CiAa^/-^ \v,2 ■ e I 
5 

66 



1/5 



+ ^C,W/^|z/r2-er + .... (7) 

The dependence of x on the impact velocity does not 
allow to derive from the Boltzmann equation a time- 
independent equation for the scaled distribution function 
/. Contrary to the case of £ = const. 0, [2^ 113 the 
scaled distribution function depends explicitly on time. 
Therefore, we write for a gas of viscoelastic particles 



f{v,t) 



:mt) 



(8) 



with the number density of the granular gas n and the 
thermal velocity vxit) defined by the granular tempera- 
ture: 



dv—f{v,t) = -n^^. (9) 



Hence, the shape of the velocity distribution function, 
characterized by the rescaled function /(c, i), does not 
persist but evolves along with temperature 113 ■ We 
wish to stress that the time dependence of f{c,t) is 
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caused by the dependence of the factor x on. the im- 
pact velocity. Contrary, for a gas of simpHfied particles 
(e = const.) we obtain x = l/s'^ = const, and, therefore, 
the rescaled distribution function is time independent, 
/(c,i) = 

For slightly dissipative particles the velocity distribu- 
tion function is close to the Maxwell distribution. It 
m ay b e described by a Sonine polynomial expansion 
dElllllll: 



/(c,t) = (/.(c) l-f^ap(t)^p(c2) 



(10) 



where 0(c) = vr exp (— c^) is the scaled Maxwell dis- 
tribution, Sp{x) are the Sonine polynomials. 



So{x) - 1 
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(11) 



X" 5x 15 



etc. 



and ak{t) are the time- dependent Sonine coefficients, 
which characterize the form of the velocity distribution 
[23, 13- The first Sonine coefficient is trivial, ai = 0, 
due to the definition of temperature [2^ , while the 
other coefficients quantify deviations of the moments Icf'^ 
of the velocity distribution from the moments of the 
Maxwell distribution (0*^)0, e.g. 



02 = 



(c^>o 



(12) 



Thus, the first nontrivial Sonine coefficient 02 character- 
izes the fourth moment of the distribution function. 

For small enough inelasticity {e > 0.6) the distribu- 
tion function is well approximated by the second Sonine 
coefficient 02 p3. 123L l28j . i.e. higher coefficients ak — 
for /c > 3 may be neglected. With this approximation 
the evolution a granular gas of viscoelastic particles in 
the homogeneous cooling state is described by a set of 
coupled equations for the granular temperature and for 
the second Sonine coefficient US, ^3 : 



dT 

m 

da2 



5^2(1 + 02) 
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(13) 
(14) 



The coefficient B = B{t) = VT{t)<J^ g2{<j)n is propor- 
tional to the mean collision frequency. The moments of 
the collision integral, ^2 and /i4 read with the approxi- 
mation / = 0(c) [1 + 02(^)^2 (c^)] : 

/Xp = -i j dci j dc2 J deQ{-ci2 ■ e)\ci2 ■ e\(j){ci) 

X (C2) {1 + a2 [S2 (cl) + S2 (cl)] + al S2 (c?) ^2 (c^) } 



X A (cf 



(15) 



where 



AV- (c- ) = VXcD - ^(q) 



(16) 



denotes the change of some function -0 (q) according to 
a collision. The coefficients /ip depend on time via the 
time dependence of a2- For small enough dissipation the 
moments ^2 and /X4 may be obtained as expansions in 
the time-dependent dissipative parameter 5' , 



5'{t) = Aa^'^ [2T{t)f'^^ = 5 



2T(t) 



Tn 



1/10 



(17) 



with S = Ao^^^Tq^^^ and with the initial temperature 



To . These expansions read |2J, |23| 

2 2 

k=a n=0 

2 2 



(18) 



k=0 n=Q 



where Akn and Bkn are numerical coefficients. They may 
be written in the compact matrix notation (rows refer to 
the first index): 
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(19) 



12500^0 
4459833 



70400000 



UJi j 



with 



Wo = 2\/2^2i/i"r ( « 6.485 



V2^2i/5r cl « 9.285 . 



(20) 



The coupled equations H13|) and H14|) together with H17|l . 
(|18|l and H19|l determine the evolution of a granular gas 
of viscoelastic particles in the homogeneous cooling state. 
In particular, they define the velocity distribution func- 
tion which is the starting point for the investigation of 
inhomogeneous gases. 



B. Adiabatic approximation for the second Sonine 
coefficient 



In the limit of small dissipation, 5 <C 1, the coupled 
equations (|13ll . (|14|) may be solved analytically. In linear 
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app roximation with respect to S the solution reads |2i 



To V "^0 



'5/3 



where we introduce the characteristic time 



-qo6Tc{Oy 



(21) 



(22) 



with the initial mean collision time Tc(0) and the constant 



90 



2i/5r ^« 0.173, 

\10J 8 



(23) 



Correspondingly in linear approximation the second So- 
nine coefficient depends on time as 24] 

a,{t) = -^w{tr^ {Li [w{t)] - Li [wiO)]} , (24) 
5 



with 



w{t) = exp 



(90^)"' ( 1 



To 



1/6- 



and with the logarithmic integral 



Li(a;) 



1 



dt. 



(25) 



(26) 



For small dissipation i5 the time dependence of 02 re- 
veals two different regimes: (i) fast initial relaxation on 
the mean collision-time scale ~ 7"c(0) and (ii) subsequent 
slow evolution on the time scale ^ tq':^ '''c(O), i.e. on the 
time scale of the temperature evolution. Therefore the 
coefficient a2 (and hence the form of the velocity distri- 
bution function) evolves in accordance with temperature. 

For the hydrodynamic description of granular gases we 
assume that there exist well separated time and length 
scales. The short time and length scales are given by 
the mean collision time and the mean free path, and the 
long time and length scales are characterized by the evo- 
lution of the hydrodynamic fields (to be defined in the 
next section) and their spatial inhomogeneities. The hy- 
drodynamic approach corresponds to the coarse-grained 
description of the system, where all processes which take 
place on the short time and length scales are neglected. 
Therefore the first stage of the relaxation of the velocity 
distribution function does not affect the hydrodynamic 
description and only the second stage of its evolution on 
the time scale tq » Tc(0) is to be taken into account. 
To this end we apply an adiabatic approximation: We 
omit the term da2/dt in the left-hand side of Eq. H14|) 
which describes the fast relaxation and assume that a2 is 
determined by the current values of /Z2 and 1^4 due to the 
present temperature. Hence, in the adiabatic approxima- 
tion 02 is determined by 



5/i2 (1 + a2) - /i4 = . 



(27) 



Using (|18|l for ^2, find 02 as an expansion in the 

small parameter S': 

a2 = a2iS' + a22S'^ + ... (28) 
3ujo 



021 = 



a22 



-0.388 



20V27r 

12063 cjg 

640000 V ^ 40 



27 LUi 



2.752. 



Figures n and 121 show 02 in adiabatic approximation due 
to Eq. (|28|l together with the numerical solution of Eqs. 
(|13I14|I and the result of the hnear theory, Eas. (|24l25|l . 

The adiabatic approximation is rather accurate for small 
dissipation, 6 < 0.05, after the initial relaxation (5-10 
collisions per particle) has passed. Even for the larger 
dissipation {S = 0.2) it is in agreement with the numer- 
ical result. This value of the dissipative parameter cor- 
responds to the initial coefficient of restitution e « 0.75 
for the thermal velocity. The adiabatic approximation 
becomes more and more accurate as the system evolves. 
Contrary, the linear theory fails with increasing S and be- 
comes qualitatively incorrect for larger 5. In Figs, ^and 
|21the evolution of the system has started from a Maxwell 
distribution, i.e., 02(6) = 0. Fig. |3| shows the relaxation 
of 02 from certain initial conditions 02 (0) 7^ which cor- 
respond to non-Maxwellian distribution functions. The 
relaxation occurs during the first 5-10 collisions per par- 
ticle. Then the adiabatic approximation becomes valid 
even for intermediate values of the dissipative parameter 
6. The adiabatic approximation yields very accurate re- 
sults for the evolution of temperature. It coincides almost 
perfectly with the numerical solution for the short and for 
the long time scales and for all values of the dissipative 
parameter < i5 < 0.2 (see Fig. 0)). Thus we conclude 
that the adiabatic approximation may be applied for the 
hydrodynamic description of granular gases. By means 
of Eq. for the second Sonine coefficient we can de- 
termine the moments /X2, /X4, and the cooling coefficient 
C: 



2 2 2T 



where 



3 V m 



9 2 /2 



(29) 



(30) 



Below we will need the derivatives of 02 and C with 
respect to temperature and density. Using 

^dS' _ S' 

~ To 

according to the definition of of (5', Eq. ()17|l . we obtain 



T- 



(31) 
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(32) 
(33) 
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FIG. 1: Evolution of the second Sonine coefficient 02 in the 
homogeneous cooling state on the short time scale. Solid line 
- numerical solution of Eas. ll3ll4L dashed line - adiabatic 
approximation 12811 . dot-dashed line - the result of the linear 
theory, Egs. 1241251 . The time is given in collision units Tc(0). 
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FIG. 2: The same as in Fig. Q but for longer time. The 
adiabatic approximation improves as the system evolves. 



III. HYDRODYNAMIC EQUATIONS AND 
TRANSPORT COEFFICIENTS 



A. Derivation of hydrodynamic equations 



For the derivation of hydrodynamic equations from the 
Boltzmann equation it is assumed that there exist well 
separated time and length scales. As already briefly dis- 




FIG. 3: The same as in Fig. Q but for different initial values 
of a2. The fast relaxation takes place during the first 5-10 
collisions, independently on 02(0). 



FIG. 4: Evolution of temperature in the homogeneous cooling 
state on the short time scale (top) and for longer time (bot- 
tom). The notations are the same as in Fig. [H The dotted 
line shows the asymptotic power-law &t t ^ 00. 



cussed in the preceding Section, in granular gases there 
are (at least) two sets of scales: The microscopic scales 
are characterized by the mean collision time and the 
mean free path. The macroscopic scales are given by 
the characteristic time of the evolution of the hydrody- 
namic fields and the size of their spatial inhomogeneities. 
Hence, scale separation means that the macroscopic fields 
vary very slowly in space and time if measured in the mi- 
croscopic units. This assumption allows for a gradient 
expansion in space and time, i.e. the application of the 
Chapman-Enskog method ^ . 

The application of the Chapman-Enskog approach to 
granular gases is more sophisticated than its application 
to molecular gases: For molecular gases the unperturbed 
solution (the basic solution) for the velocity distribution 
function is the Maxwell distribution while for granular 
gases the basic solution is the time-dependent velocity 
distribution of the homogeneous cooling state. 

For gases of particles in the homogeneous cooling state, 
which collide with e = const., the velocity distribu- 
tion function depends on time only through the time- 
dependent characteristic velocity. The shape of the ve- 
locity distribution function, given by /, is fixed by e 



and does not depend on time. Hence the evolution of 
temperature describes the evolution of the velocity dis- 
tribution function exhaustively. For granular gases of 
viscoelastic particles the shape of the velocity distribu- 
tion function, characterized by the Sonine coefficients a^, 
is time dependent due to the time dependence of these 
coefficients. In the previous section we exemplified this 
property for the first non-trivial coefficient 02 (t) , assum- 
ing that for small dissipation higher coefficients may be 
disregarded. Therefore, the evolution of granular gases 
of viscoelastic particles is described by the time depen- 
dence of temperature, i.e., by the second moment of the 
velocity distribution function, and by the time depen- 
dence of its higher-order even moments which are char- 
acterized by the Sonine coefficients 02, 03, etc. [see Eq. 
(|12|l ]. Hence, the hydrodynamic description of granu- 
lar gases whose particles collide with e — e(g) requires 
an extended set of hydrodynamic fields which includes 
the higher-order moments. A regular approach to derive 
hydrodynamic equations for this extended set of fields 
is the Grad method |J|. Alternatively, for small dissi- 
pation an adiabatic approximation can be applied: It 
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is assumed that the shape of the velocity distribution 
function, albeit varying in time, follows adiabatically the 
current temperature. Correspondingly, the higher-order 
moments are determined by the temperature too. With 
this approximation a closed set of hydrodynamic equa- 
tions for density n{f, t), velocity u{f, t), and temperature 
T(r, t) may be derived. These fields are defined respec- 
tively by the zeroth, first, and second moments of the 
velocity distribution function: 



(35) 



n (r, ^) = y ^'^Z ^' 

n (r, t) u (r, t) — J '^''^ ^/ ('^i i 

3 f 1 

-n (r, t) T (f, t)= dv -mV'^ f (r, v, t) 



where V = v — u{r,t). For small dissipation the veloc- 
ity distribution is well approximated by the second So- 
nine coefficients 02 (t), i.e., higher-order coefficients are 
neglected. In adiabatic approximation a2 is determined 
by temperature according to (|28|l . 

For simplicity of the notation we consider a dilute gran- 
ular gas and approximate the Enskog factor 52(0') ~ 1- 
Multiplying the Boltzmann equation for an inhomoge- 
neous gas 

d_ 

dt 



+ i/i-V /(f, 1/1,0 = J (/,/), (36) 



correspondingly by fj, wi, rnui/2, and integrating over 
dvi we obtain the hydrodynamic equations (see e.g. 

Qti -> 

— +V-(™i) = 0, (37) 
+ M- Vu+(nTO)"iV-P = 0, (38) 



du 
'dt 



u-VT + — (P -.Vu + V ■ q) +CT = 0. (39) 



dt 



3n 



The cooling coefficient C in the sink term C.T may be 
written as 



UnT 



j dvi J J deO {-V12 ■ e)\vi2 ■ e\ 

X / (r, t/i, t) / (f, V2,t) [vi2 ■ ef (I - e^) . (40) 
The pressure tensor P and the heat flux q are defined by 

(r, = / (^) / (r, t/, t) di/ + pS,j (41) 



<f(f,<) = / S{V f{r,v,t)dv, 



(42) 



where p = nT is the hydrostatic pressure. The velocity 
tensor Da and the vector 5" read 



siv] . {!^-It]v. 



(43) 
(44) 



The structure of the hydrodynamic equations, except for 
the cooling term (T, coincides with those for molecular 
gases. 



B. Chapman-Enskog approach 



The system (|37I38I39II is closed by expressing the pres- 
sure tensor and the heat flux in terms of the hydrody- 
namic fields and the fields gradients. To this end we 
apply the Chapman-Enskog approach This method 
is based on two important assumptions: The evolution 
of the distribution function is completely determined by 
the evolution of its first few moments, i.e., it depends on 
space and time only through the hydrodynamic fields: 



/ {r, v,t) ^ f [v, n (f, t) ,u{r,t),T [r, t)] 



(45) 



As the second precondition it is assumed that the gas 
is only slightly inhomogeneous on the microscopic length 
scale which allows for a gradient expansion of the velocity 
distribution function: 



(1) 



X2f(2) 



(46) 



where each power k of the formal parameter A corre- 
sponds to the order k in the spatial gradient. Thus, 
refers to the homogeneous cooling state, Z^^-* corresponds 
to the linear approximation with respect to the fields gra- 
dients, f^"^^ is the solution with respect to quadratic terms 
in the field gradients, etc. With these assumptions the 
Boltzmann equation may be solved iteratively for each 
order in A, together with the hydrodynamic equations for 
the moments of the velocity distribution function. The 
solution in zeroth order in A yields the velocity distri- 
bution function for the homogeneous cooling state /'"^ 
and the corresponding evolution of temperature. The 
function is then used to compute P and (f, yielding 
P'-*''' = 6ijp^ (f'-*'-' — and the hydrodynamic equations 
for the ideal fluid. These first-order equations contain 
only linear gradient terms. Then /(^•' may be found em- 
ploying the first-order hydrodynamic equations and the 
distribution function /*^°^ The obtained f^^^ as well as 
the corresponding expressions for P*-^) and g^^^ are linear 
in the field gradients: 



Py ^p5i3 - ri ( ViUj + V-jUi - - %V • u ^ ^^^^ 



q 



The transport coefficients 77, k and fi in these equa- 
tions are expressed in terms of Z'^'. Hence, within the 
Chapman-Enskog approach for each order in the gradi- 
ent expansion a closed set of equations may be derived. 
Keeping only first-order field gradients for the distribu- 
tion function and respectively for the pressure tensor and 
for the heat flux as in Eq. (|47|l . the Navier-Stokes hydro- 
dynamics is obtained. Keeping next-order gradient terms 
corresponds to the Burnett or super-Burnett description. 



8 



We will restrict ourselves to the Navier-Stokes level |32| 
and skip for simplicity of the notation the superscript 
"(1)" for P and q in the above equations. 

The Chapman-Enskog scheme also assumes a hierar- 
chy of time scales and respectively a hierarchy of time 
derivatives: 



dt~ dt ^ dt ^ dt 



(48) 



(fc) 



where each order k in the time derivative, d^''^ /dt = d'^ 
corresponds to the related order in the space gradient. 
Consequently, the higher the order in the space gradient, 
the slower is the according time variation. Using the for- 
mal expansion parameter A we can write the Boltzmann 
equation, 



1^ + ^^ 



(/(°)+A/W + ...),(/(")+A/W + ...)] (49) 



and collect terms of the same order in A. The equation 
in zcroth order 



dt 



(50) 



coincides with Eq. © for the homogeneous cooling state. 
According to the Chapman-Enskog scheme we obtain the 
velocity distribution function in zcroth order: 



/^°^ ^ = ^^lA [1 + (C^)] 0(C) 



r,t) 



(51) 



where c — {v ~ u) /vt — V/vt- The corresponding hy- 
drodynamic equations in this order read 



dt 



■n = 0, 



a(o) 

dt 



■u = 0, 



dt 



r = -c(°)r, (52) 



where C^"-' is to be calculated using Eq. H40() with / = 
/(''^. In this way we reproduce the previous resuh if^ . 
Collecting terms of the order 0(X) we obtain. 



5(0) fii 
dt 



where we introduce 



= 0, 
(53) 

(54) 



The corresponding first-order hydrodynamic equations 
read 



9(1 



dt 

dt 

~dr 



V {nu) 
u ■ Vu - 
u-VT- 



1 

nm 
2 
3 



TV -u- 



(55) 



The term C^^^ is found by substituting / = -I- A/f^) 
into Eq. H4()|l and collecting terms of the order 0{X). 
These terms contain the factor {v, r, t) /'^^ (t/, r, t) in 
the integrand. They vanish upon integration according 
to the different symmetry of the functions /^"^^ and f^^\ 
thus, C^^^ = [23. 

Since the distribution function /("' is known, we can 
evaluate those terms in (j53|l which depend only on Z^"^: 



dn 



du \ dt 

dfio) /ad 



dt 
vi ■ Vm 



vi ■ Vn 



dT V dt 



■T + vi- VT 



(56) 



With vi ^ V + u, the time derivatives of n, u and T given 
in H55|) . and the relations 



dn n ' 



a/(o) 5/(0) 



du 



dV 



(57) 



which follow from Eq. H51|) . we recast Eq. H53|) for /(^^ 
into the form: 



a(o) f{l) 

dt 



(/(o)j(i)) =/(o) (v-u-f -Vlogn) 
a/(o) 



5T 



-TW -u-V -WT 



F • V U 



(58) 



where we take into account C'^^ = 0. 

The right-hand side is known since it contains only 
/C'). It is convenient, however to rewrite it in a form 
which shows explicitly the dependences on the fields gra- 
dients. Employing 

IT T 

V,p= — VaogT+— V.logn, (59) 

nm TO TO 



the right-hand side of H58() yields 



dt 



+ J(1) ^fiO)j(l) 



^ A-VlogT + B - Vlogn + C^jVjU^, (60) 



with 



A{V) 



-VT 



B{v) = -y/(o) 



Qf{Q) QfjO) 

dT m Qv 

TO dV 



2T a/(o) 



(61) 
(62) 
(63) 
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We will calculate these terms below. From the form of 
the right-hand side of Eq. l|t)U|l we expect the form of its 
solution 



/(^^ = a • VlogT + /3- Vlogn + 7yVj? 



(64) 



which is the most general form of a scalar function, which 
depends linearly on the vectorial gradients VT, Vn and 
on the tensorial gradients V jUi. The coefficients a, (3 and 
7jj are functions of V and of the hydrodynamic fields n, 
u and T. 

We derive now equations for the coefficients a, (3 and 
7y by substituting /^^^ as given by Eq. into the 

first-order equation (|50|l and equating the coefficients of 

the corresponding gradients. To this end we need df^ f^^'' 
and, therefore, the time derivatives of the coefficients a, 
(3, and 7^ : 

9(0)a dad'-"^T dad^^^n da d^^^T 



dt 



dT dt 



dn dt 



dui dt 



= ~C'''T^^ (65) 

where we use Eqs. (|52|l in zeroth order. Similarly we 
obtain 



dt ^ dT' 



^ = _^(o)r^ (66) 



dt ^ dT 

and, respectively, the time-derivatives of the gradients 

0(0) 

dt 



-V logn =0 



^VlogT = -VC(°) 



(67) 



The derivatives of C'"' are given by Eqs. (|33I34|I . From 
(|^ . (|^ . and (EZI) we obtain: 



where J^^^ is defined by Eq. (|54|l . 



C. Kinetic coefficients in terms of tlie velocity 
distribution function 

From the definition of the pressure tensor Eq. 141(1 
and its expression in terms of the field gradients Eq. I(47|l 
follows 



(/^°^ + a • V log T + ^ • V log n + 7fc; V,Wfe) dV 



= -f] [ V^uj + VjM, - -(5y V • w ) , (72) 



where the tensor £>, 



has been defined above. The 



integrals 

Jo^jf^'^UV^ j D.^ddV^ jD,jPdV = (73) 

vanish since Dij is a traceless tensor and J*^"^ depends 
isotropically on V. Moreover, as will be shown below, 
the vectors a and /3 are directed along V, hence, the 
respective integrands are odd functions of V. Therefore, 
only the term with the factor ^ki^iUk in the left-hand 
side of Eq. 1(721) is non-trivial. Equating the coefficients 
of the gradient factor ViUk we obtain 



D.jjkidV = -1] [ 6u6kj + SijSkt ~ -SijSki ) • (74) 



For k ^ j, I = i the last equation turns into 



Dijjj^dV = -rj ( SuSjj + SijS^j - -SijS.^j j = -IO77, 

(75) 

{SiiSjj = 9 and SijSij ~ 3 according to the summation 
convention) and yields the coefficients of viscosity 



dt 



dT 



(68) 



where we use Eo. dMll . If we insert into Eq. 

and equate the coefficients of the gradients we arrive at 
a set of equations for the coefficients a, (3 and 7ij : 



A 



B 



(69) 
(70) 



_C(o)r^ + j(i)(/(o),^,^.) ^c,, (71) 



V = ^^jD,,{v)j,,{v)dV. (76) 

Using Eq. ((41() for the heat flux and its corresponding ex- 
pression in terms of the field gradients, Eq. I(47() . we can 
perform a completely analogous calculation and arrive at 
the kinetic coefficients n and fi: 



^ = -^J dVS[v)-a[v) (77) 



with S defined by Eq. ((44(1 . The coefficient of thermal 
conductivity k has the standard interpretation, while the 
other coefficient /i does not have an analog for molecular 
gases. 
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IV. COEFFICIENT OF VISCOSITY 



where we take into account that 



The viscosity coefficient is related to the coefficient ^ij 
[see Eq. ((7^ ] which, in its turn, is the solution of Eq. 
(|71|l with the coefficient Cij in the right-hand side. Let 
us first find an explicit expression for Cij . 

According to Eq. H51|l the velocity distribution func- 
tion depends on temperature through the thermal veloc- 
ity vt and additionally through the second Sonine coeffi- 
cient a2- Hence the temperature derivative of f^^^ reads 



df 



(0) 



dT 



with /a/(V) being the Maxwell distribution 



/m = -^<l>{c) with cf){c) 



1 



r3/2 



exp {—(? 



Using the relation 



dVi V dV 
and Eq. (|79|l the coefficient Cij reads 

1 „ 1 



(79) 



(80) 



(81) 



C^J [ V 



With 

1 0/(0) _ 
V dV ^ 

we obtain 

Cij 



ViVi - -SiiV^ , 
y^yJ 3"*^" J V dV 



(82) 



[1 + a,S, (c^)] /m + ^(c2-0 a,fM 

(83) 



1 

— A, 

rp '■J 



l + a2\S2 ic 



hi{v) 



+ ls,,S2{c')fM{V)T^, (84) 



where Dij (j^^ has been defined in Eq. (|43|l . 

The expression for Cij determines the right-hand side 
of Eq. H71|) for 7^ and hence it suggests the form for 
7y . For small dissipation (when a2 is small) and for 
small fields gradients we keep only the leading terms with 
respect to these variables. Therefore, we seek for 7^ in 
the form 

l^,{v)=YDu{v)fM{V), (85) 

where 70 is a velocity- independent coefficient, i.e., we 
neglect the dependence of jij on 02 [s^. The viscosity 
coefficient Eq. (|76|l reads then 

/ '^^^vDijIm ^ -lonT , (86) 



(87) 



according the definition of Dij , Eq. 143(1 , with the sum- 
mation convention. Moreover, we have used the fourth 
moment of the Maxwell distribution. From Eq. (|86|l fol- 
lows 



70 



JL 

nT 



Multiplying Eq. H71|l by Dij{Vi), integrating over Vi 
and using Eq. l(7^ yields 

+ J dV^iA, (l?i) j('^(/^°\7y ) • (89) 

To evaluate the first term in the right-hand side we use 
Eq. (|63|l for Cij, the relation 

^A. ^ (V.V, - = mV, (1 + 1^.,) , 

(90) 

the definition of temperature, Eq. H35() . and notice that 
DijSij = since Dij is a traceless tensor [see Eq. ()43(l ]. 
Integration by parts then yields 



8 2T f r^f^") 



(91) 
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lOnT. 



For the second term in the right-hand side of H89|l we use 
the definition of J^^\ (1541). and obtain 



= - / dyiA,/(/(°\7., ) - / dV?iA, /(7.„/^°^) 



(92) 



We apply the property of the collision integral 

j dl^iA,/(/("\7.,) - j dV,D,,l(j,,J^°^) 
= Y y M y" dV2f^"^ [Vi) 7y (^2) J dee {~Vi2 ■ ' 



X Vi2-e A D,, [VA+ A, [V2) , (93) 
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where Ai/; (vi) = ^ (v/) — ^(vi) denotes as previously the 
change of some quantity tpi^i) due to a coUision. Equa- 
tion l|92|l turns then into 



dViD, 



1^J {V2) j dee (-^12 • e ) 



V12 ■ e 



D, 



(Vi) + (^2)] . (94) 



We write the factors in the last integral using the dimen- 
sionless velocities V1/2 — vtCi/2- 



^mv^Dij (c) = mv^ (ciCj - 



T \v. 



(95) 



3 I mv^D^j (c) 4> (c) 



7., (v) 

and recast Eq. (|94|) into the form 

dViD.j (^i) J(i) (/(°\ 7^,) = 4.T]VTna^n^ , (96) 

where we substitute 70 — —rj/nT and where fl^ is a nu- 
merical coefficient defined by 

^l^ = j dci J dc2 j dm (-C12 • e) |ci2 • e I /(°) (ci) 

X (C2 ) (C2 ) A (ci ) -f (C2 )] . (97) 

The coefficient $7^ may be expressed in terms of the dis- 
sipation parameter 5' and the second Sonine coefficient 
02 (see Appendix): 



with 



97 679 

ui^ =^^^1 I 02 

' 165 44000 



(98) 



(99) 



Substituting Eqs. H91I96I98I) into Eq. ^ and using Eq. 
(El for C^"^ 

we obtain an equation for the coefficient of 

viscosity rj: 



{.,6'-.2nT^ 



= ^-{nj, + w,5' -W2nv-~2^- (100) 
We seek the solution as an expansion in terms of 5': 

?7 = 77o(l + (5'r)i+<5'2r)2 + ...) . (101) 



The solution in zeroth order 



16ct2 



mT 



(102) 



is the viscosity coefficient for a gas of elastic particles 
(Enskog viscosity), while the coefficients fji and 772 ac- 
count for the dissipative properties of viscoelastic parti- 
cles. With Eq. (|31|l the temperature derivative of the 
viscosity coefficient reads 



T 



drj 
df 



(103) 



We substitute Eqs. (|101I103|I into Eq. pOO|) . express 02 
in terms of 6' according to Eq. H28I) and collect terms of 
the same order in 6' . This yields the equations for fji, 772, 
etc., whose solutions read 



: 0.483 



359 V2^ 
=3840^'^" 

41881 _ 567 c.iV2^ _ 
^2304000 TT 28160 tt ~ ' 



(104) 



with Wq/i given by Eq. H20|) . 

Thus, we arrive at the final expression for the viscosity 
coefficient for a granular gas of viscoelastic particles: 



-^J—(1 + 0A83 5' + 0.0MS'^ + ...) . (105) 

16(7^ V TT 



In contrast to granular gases of simplified particles (e = 
const.), where ij cx VT, for a gas of viscoelastic particles 
there is an additional temperature dependence due to the 
time-dependent coefficient 5'. 



V. COEFFICIENT OF THERMAL 
CONDUCTIVITY AND THE COEFFICIENT fi 

To find coefficients k and /i we need the coefficients 
a and /3 which are the solutions of Eqs. (|69I70|I . The 
functions A and B in the right-hand sides may be found 
from Eqs. H61I62() : 



A 



-^S{v) [l + a2 (52 (c2) + 1 - c^)] /m 



VS2 (c2 



0-21 Q22 

To" T 



fM (106) 



(107) 



Keeping only leading terms with respect to the gradients 
and 02, we choose a in the form 



d = -^s{v)fM{v) 



(108) 



12 



with ai being the velocity independent coefRcient. This 
Ansatz for a yields the coefBcicnt of thermal conductiv- 
ity, 



1 

'3T 



dVS (y) • a (f) 
dV 



ai 

3r2 



which implies 



ai 



2m 
f 

5nT 



5 nT 
= (109) 
2 m 



(110) 



Multiplying Eq. ^ for a by S* [Vi^ j T, integrating 
over Vi and using Eq. (|77|l for k we obtain 



T 



T 



(111) 

To evaluate the first term in the right-hand side we use 
Eq. (gH) for S {V^ and Eq. (dH) for A {v^ : 



1 

T 



rives' 



(^).i(v) 



■2 / f c2 



- - ) c2 [1 + a2 (^2 (c^) + 1 - c2)] /m 



^2 (c^) /m . 















L 10 


5 



(112) 

With the Maxwell distribution Eq. (|80|) the first term in 
the last equation H112(l reads 



+a2 



cH{c)\c'-- \ dc 



cU(c) 





15 



dc 



-^"2 15 2 15 n-T 

^ —nvrp + —nvrpa2 = — (H-2a2), (113) 

4 2 2 m 

where integration over the angles has been performed and 
the integral 



r ( 2^ 2fcj (2fc-l)!! _ 
exp {-x^) x^'^dx = ^ 



(114) 



was used. Very similar calculations give the second term 
in Eq. ifTHI) : 



15 



V 10 



Q22 



(115) 



Summing up Eqs. (|113|l and (|115|l we obtain the first 
term in the right-hand side of (|lll(l : 

1 f jT-> o 15 «r / 21 ^, 11 ^,2\ 

(116) 

The second term in the right-hand side of Eq. 1)111(1 may 
be again written using the basic property of the collision 
integral [see Eq. (|^ ]: 



1^12 ■ e 



Using the dimensionless variables 



(117) 



a \ V) = -ai VT\-^] 4'{c)S (c ) 



(118) 



and Eq. 1(110(1 for ai, we recast the last equation into the 
form 



1 
T 



dVs(yyj^^'^ =--nvTna'^n^. (119) 

The coefficient D,^ is defined by 

17, = y" dci y" dc2 j deO (-C12 • e) |ci2 • e1 /(°) (ci) 

X 0(C2)5(C2) • A [^(ci)-|-S'(c2)] . (120) 
This coefficient reads (see Appendix) 



(121) 



where 



uo =4^2^ (-^ + 3^ 

/17 9 
U-500"^ 
/1817 1113 
-"'^ 1,140" "352000"^ 



(122) 



Substituting Eqs. ((116I119I121() into Eq. lfTTT|l and using 
Eq. ((^ for C'-^-' we arrive at 

T^kT^/^ {^,5'-u25'^ + ...) 

^-kT^'^ (uo + 5'm-5'^U2 + ...) 
5 

15T3/2 nr / 21 ^, 11 ^,2^ 



(123) 
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We solve this equation with the Ansatz 

K = KO (1 +(5'ki + (5'2k2 + . . .) 

where 



Kq 



75 
64cr2 



(124) 



(125) 



is the Enskog thermal conductivity for a gas of elastic 
particles. Substituting Eq. 11241) into Eq. H123|l and 
equating terms of the same order in 6' we obtain the 
coefficients 



Ki 



K2 



487 VSttwo ^ 
6400 n 
1 / 2872113 



0.393 

, 78939 



TT V 51200000 



OJ, 







140800 



2'KiLj'l 



4.904. 

(126) 



Hence, the coefficient of thermal conductivity for a gran- 
ular gas of viscoelastic particles reads in adiabatic ap- 
proximation 



4.904(5' 



(127) 



Similar to the viscosity coefficient, the coefficient of ther- 
mal conductivity of a granular gas of viscoelastic particles 
reveals an additional temperature dependence as com- 
pared with gases of simplified particles where e — const. 

The evaluation of the coefficient fi may be performed 
in the same way as k, i.e. choosing /3 in the form 



0=-Y^{v)fM{v). 



(128) 



with the velocity independent coefficient p. The only 
difference is that the expansion of /i in terms of the dissi- 
pative parameter S' lacks the term in zeroth order since 
fi vanishes in the elastic limit. Since the calculations are 
completely analogous to that for n, we present here only 
the final result: 



with 



19 cjoV27r 
Mr « 1.229 

, 1 / 58813 2 1 \ 1 .1. 

112 =- ojn v2ttuji « 1.415 . 

^ TT V 640000 " 40 / 



(129) 



(130) 



KqT 



1.229 5' + 1.415 5'2 + ...) 



Finally, using the coefficients 
2m 

b^^' 5T 



Ql 



2m 



70 



nT 



(131) 



ry, (132) 



(the result for /3i may be derived analogously as for ai) 
in the relations I|108I128I85|I for a, (3, jij we obtain an 
expression for the first-order distribution function /'^^ 
which depends linearly on the field gradients according 
to Eq. 

VI. TWO-DIMENSIONAL GRANULAR GAS 

So far we have restricted ourselves to three-dimensional 
systems although the calculations are identical for gen- 
eral dimension d. Of particular interest is the case d — 2. 
Since molecular dynamics simulations are frequently per- 
formed for two-dimensional systems, in this section we 
present the results for two-dimensional gases of viscoelas- 
tic particles. We wish to stress that these systems are, in 
fact, quasi two-dimensional, since we still use the coeffi- 
cient of restitution Eq. for colliding spheres. Hence, 
we assume that the motion of the spherical particles of 
the gas is restricted to a two-dimensional surface. 

The hydrodynamic equations for two-dimensional 
gases have the form 



dn 
'dt 



+ V • {nu) = 



^-f u- W+(nm)-^V-P = (133) 
ot 

^ + u • VT + i (p,,v,u, + V • (f) + cr = , 

with the pressure tensor and the heat flux 

Pij = nTSij - 77 (\/^Uj + VjU, - S,jV ■ uj 
q = —kVT — fiVn . 
Correspondingly, the transport coefficients read 

7? = ^70 (1 + S'fii +S'^fi2-\ ) , 

where 6' (t) ^ 5 [2T (<) /Tof^^" and 



(134) 



(135) 



riQ 



1 

2^ 

29 



mT 

TT 
2^ 



V2 

Similarly, 



640 TT 
111 



-U!o 



0.234 

569 ujiV2^ 



(136) 



160000 TT 14080 



0.308 . 



Kq (1 + S'ki + S''^K2 + • • • ) 



(137) 



Thus, the coefficient fi reads in adiabatic approximation with 



Kq 



-V- 

a V TTUi 
433 V2ttujo 



3200 TT 
1 ( 1749573 



= - 



TT \ 12800000 



-0.700 

95619 
70400 



2ttijJ^ 



(138) 
11.89. 
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and 



^(<5'Ai + <5'^A2 + ---) 



Ml 



n 

7 wo\/27r 
20 

1 / 7411 



1.811 



(139) 



80000 



Wo + ^\/2^wi^ w 0.056. 



The cooling rate C*-"^ reads for a two-dimensional gas of 
viscoelastic particles 

V m \2 



UJ2 



1 9V27r , 
-wi + ^TTT— t^n ~ 5.246. 



(140) 



2 5007r 

The numerical constants ujq/i are given by Eq. (|20|l . 

VII. SUMMARY 

We have derived the hydrodynamics of granular gases 
of viscoelastic particles. Collisions of viscoelastic parti- 
cles are characterized by an impact velocity dependent 
coefficient of restitution. We have used the Chapman- 
Enskog approach together with an adiabatic approxima- 
tion for the velocity distribution function, which assumes 
that the shape of the velocity distribution function fol- 
lows adiabatically the decaying temperature. We have 
compared the numerical solutions for temperature T{t) 
and for the second Sonine coefficient 02 {t) with the corre- 
sponding adiabatic approximations and have found good 
agreement up to intermediate dissipation. To derive the 
hydrodynamic equations and transport coefficients from 
the Boltzmann equation we have modified the standard 
scheme to account for the time dependence of the basic 
solution. This modification takes into account the time 
dependence not only due to the thermal velocity, as for 
e = const., but also due to the evolution of the shape of 
the distribution function as given for gases of viscoelastic 
particles. 

Transport coefficients and the cooling coefficient for 
dilute granular gases of viscoelastic particles read 

r, = r^o{l + S'm+6'% + ...) (141) 

K = Ko{l + d'kl+S'^k2 + ■■■) (142) 

(143) 

nT 



^ = [S fli+S fl2 

n 



(144) 



where 770 and kq are the Enskog values for the viscosity 
and the coefficient of thermal conductivity and 5' is the 
time dependent dissipative parameter. The numerical 
coefficients ryi/27 '«i/2; Ai/2; and C1/2 are given in Table 



TABLE I: Numerical coefficients for Eqs. II141I1441 





3d 


2d 




3rf 


2d 




0.483 


0.234 


'72 


0.094 


0.309) 




0.393 


-0.700 


K2 


4.904 


11.893 


Ml 


1.229 


1.811 


M2 


1.415 


0.056 


Ci 


1.078 


1.294 


C2 


-1.644 


-2.093 



The dependence on temperature and, therefore, on 
time of the kinetic coefficients Eqs. (|141|l - (|144f) differs 
significantly from the time dependence of the correspond- 
ing coefficients for granular gases of particles which col- 
lide with a simplified coUisional model e — const. This 
leads to the serious consequences for the global behavior 
of force- free granular gases [2^ . 

Under mild preconditions the presented formalism for 
the derivation of the hydrodynamic equations and the 
transport coefficients may be applied also to gases of par- 
ticles whose collision is described by a different impact 
velocity dependence than given for viscoelastic particles. 



APPENDIX A: DERIVATION OF THE 
COEFFICIENTS AND Q.^ 

For the evaluation of the numerical coefhcient Eq. 
defined by 



VI,, ^ J dcij dc2 J dSe (-C12 • e ) \ci2 ■ e I /'"^ (ci) 

X (C2) A, (C2) A [Aj (?l) + {C2)] . (Al) 

we need the factor 

Aj(c2)A[Aj(ci)+ Aj(c2)] 

= {c[ ■ C2) + (C2 • C2) - (Ci • C2f - (C2 • C2)^ 

- (c;^ + 4^ -cj- 4) (A2) 
Similarly, for the coefficient 

17, = y dci J dc2 J deQ (-C12 • e ) |ci2 • e I (ci) 

xcj){c2)S{c2)-A[s{ci) + S{c2)], (A3) 
given by Eq. (|120|l we need 



5(c2)-A S{c,) + S (02) 



= C9 - - 



2 

2 ^ X 2 ,^ ^ ^ 2 



ici-C2)ic[f 



+ (0*2 • C2) (C2) - (Ci • C2) q - (C2 • C2) C2 



(A4) 
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With (IA2II and (IA4I) we write for the coefRcient ri„ 



= j dci J dc2 j deQ (-C12 • e ) |ci2 • e?| /("^ci) 

X (C2) (c{ • C2)^ + (C2 • C2)^ - (Cl • €2)^ 
^2 1 



^ \2 2 I 12 , II 2 2\ 

(C2-C2) - (Ci +C2 -Ci-Cgj 



(A5) 



and for the coefficient 17^ 



(ci) 



X 02 (C2) (^C^ - - j [(C{ . C2) K)' + (r- . C2) (4)^ 

- (ci • C2) - (c2 • C2) C2] . (A6) 

The pre-colhsion vefocities ci , C2 as weh as after-colhsion 
velocities c{ , C2 can be expressed in terms of the center of 
mass velocity C = (ci + C2) /2 and the relative velocity 
C12 = Cl — C2 before the coUision: 



ci=C + -C12 

C2=C - ici2 



(A7) 



c{ =(7+ici2- i(l + £) (ci2-e)e 

C2 = C- ^ci2 + i (1 + e) (ci2 • e ) e . 

The coefficient of restitution is expressed in terms of the 
relative velocity C12 by 

6^1- CiS' (t) |ci2 • e\'/' + ^Cf 5'^ [t) |ci2 • e 1'/^ + . . . 

(A8) 

The second Sonine polynomial in the distribution fimc- 
tion 



/(")(ci) = </.(ci)[l + a252(c?)] 
reads in terms of C and C12 



(A9) 
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Jc^C?2 + icL {C-C,2)~\C^ 



If we replace all factors in the integrands of Eqs. HA5IA6|) 
by the corresponding expressions in terms of C and C12 
we observe that fi,, and U,^ may be written as a sum of 
integrals of the structure 



Jk,l,m,n,pM = J dci2 J '^^ J '^^Q i-Ci2 ■ c) \ci2 ■ 6 |^+" 

C ■ C12J [C ■ e j (ci2 • ef . 

(All) 

The solution of this integral in general dimension d reads 
m for n = 

Jk,Mp.a = (-l)P [1 + (-1)"] 2'+"+P+"+lf7-l 

X Pp+a+lf^mJk+mJl+m+p+a+1 , (Al2) 

for n = 1 

JkMp^c = (-1)^+' [1 + (-1)'"+^] 2'+'W"+il7-i 

X /3p+Q+2/3m+l 7fc+m+l7i+m+p+a+l ) (^13) 

and for n = 2, 

Jfe,i,™,2.p,a = {-If [1 + (-1)"] 2'+"+P+"+i 

x[(d— l)rirf] ^ 7fe+m+2 7/+m+p+Q + l 
X [(rf/3p+a+3 ~ Pp+a + l) Pm+2 + (/3p+Q+l ^ /^p+a+s) A) 



where 



0„ 



27r^ 



(AM) 



(A15) 



is the surface of a c?-dimensional unit sphere. The coeffi- 
cients 13m, 7m read 



/3m = TT 



p m+1 \ 
^ ^(<i-l)/2 ^12; 

/ m+d 'i 
- 2-™/2 V 2 ^ 



(A16) 



Following this procedure we obtain the desired coeffi- 
cients as given by Eqs. (|98() and H121() . 

The evaluation of the sums which lead to the factors fi,, 
and U,^ is straightforward, however, very lengthy. They 
can be calculated by symbolic algebra 



[1] D. Kolymbas, ed., Constitutive modelling of granular ma- 
terials (Springer, Berlin, 2000). 

[2] P. A. Vermeer, S. Diebels, W. Ehlers, H. Herrmann, 
S. Luding, and E. Ramm, eds.. Continuous and dis- 



continuous modelling of cohesive-frictional materials 
(Springer, Berlin, 2000). 
[3] T. Poschel and S. Luding, eds.. Granular Gases, vol. 564 
of Lecture Notes in Physics (Springer, Berlin, 2001). 



16 



[4] H. Grad, Comm. Pure and Appl. Math. 2, 331 (1949). 

[5] S. Chapman and T. G. Cowling, The Mathematical The- 
ory of Non-uniform Gases (Cambridge University Press, 
Cambridge, 1970). 

[6] C. K. K. Lun, S. B. Savage, D. J. Jeffrey, and N. Cliep- 
urniy, J. Fluid Mech. 140, 223 (1984). 

[7] J. T. Jenkins and M. W. Richman, Archives for Particle 
Mechanics and Analysis 87, 355 (1985). 

[8] A. Goldshtein and M. Shapiro, J. Fluid Mech. 282, 75 
(1995). 

[9] N. Sela and I. Goldhirsch, J. Fluid Mech. 361, 41 (1998). 
[10] J. J. Brey, J. W. Dufty, C. S. Kim, and A. Santos, Phys. 

Rev. E 58, 4638 (1998). 
[11] J. J. Brey, M. J. Ruiz-Montero, and D. Cubero, Phys. 

Rev. E 60, 3150 (1999). 
[12] V. Garzo and J. W. Dufty, PRE 59, 5895 (1999). 
[13] J. J. Brey and D. Cubero, in 0|, p. 59. 
[14] W. Goldsmit, Impact: The theory and physical behavior 

of colliding solids (Arnold, London, 1960). 
[15] F. G. Bridges, A. Hatzes, and D. N. C. Lin, Nature 309, 

333 (1984). 

[16] G. Kuwabara and K. Kono, Jpn. J. Appl. Phys. 26, 1230 
(1987). 

[17] R. Ramirez, T. Poschel, N. V. Brilliantov, and T. Schwa- 

ger, Phys. Rev. E 60, 4465 (1999). 
[18] H. Hertz, J. f. reine u. angewandte Math. 92, 156 (1882). 
[19] N. V. Brilhantov, F. Spahn, J.-M. Hertzsch, and 

T. Poschel, Phys. Rev. E 53, 5382 (1996). 



[20] T. Schwager and T. Poschel, Phys. Rev. E 57, 650 (1998). 
[21] S. E. Esipov and T. Poschel, J. Stat. Phys. 86, 1385 
(1997). 

[22] T. P. C. van Noije and M. H. Ernst, Granular Matter 1, 
57 (1998). 

[23] N. V. Brilliantov and T. Poschel, Phys. Rev. E 61, 2809 
(2000). 

[24] N. V. Brilliantov and T. Poschel, Phys. Rev. E 61, 5573 
(2000). 

[25] R. Ramirez, D. Risso, R. Soto, and C. P., Phys. Rev. E 

62, 2521 (2000). 
[26] C. K. K. Lun and S. B. Savage, Acta Mechanica 63, 15 

(1986). 

[27] N. Brilhantov and T. Poschel, in p. 100. 
[28] M. Huthmann, J. Orza, and R. Brito, Granular Matter 
2, 189 (2000). 

[29] N. V. Brilliantov and T. Poschel, Phil. Trans. R. Soc. 

Loud. A 360, 415 (2001). 
[30] T. Poschel and N. V. Brilliantov, preprint (2003). 
[31] I. Goldhirsch, in 0, p. 79. 

[32] It has been shown that some processes in granular gases 
require the Burnett level to derive a consistent set of 
equations 0,1^. 
[33] This may be directly checked using Eq. 1641 for f^^^ 
[34] In a more general approach yij is represented as an ex- 
pansion in orthogonal polynomials, with 1851 being the 
first order term laj. 



